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Abstract. We propose a novel approach to continuum modelling of dynam- 
ics of crystal surfaces. Our model follows the evolution of an ensemble of 
step configurations, which are consistent with the macroscopic surface pro- 
file. Contrary to the usual approach where the continuum limit is achieved 
when typical surface features consist of many steps, our continuum limit is 
approached when the number of step configurations of the ensemble is very 
large. The model is capable of handling singular surface structures such as 
corners and facets and has a clear computational advantage over discrete 
models. 

Keywords, continuum modelling, multi scale modelling, step flow, surface evo- 
lution. 



1. Introduction 

The behavior of classical physical systems is typically described in terms of equa- 
tions of motion for discrete microscopic objects (e.g. atoms). The dynamics of the 
microscopic objects is usually very erratic and complex. Nevertheless, in many 
cases a smooth behavior emerges when the system is observed on macroscopic 
length and time scales (e.g. in fluid flow through a pipe). A fundamental problem 
in physics is to understand the emergence of the smooth macroscopic behavior of 
a system starting from its microscopic description. A useful way to address this 
problem is to construct a continuum, coarse-grained model, which treats the dy- 
namics of the macroscopic, smoothly varying, degrees of freedom rather than the 
microscopic ones. The derivation of continuum models from the microscopic dy- 
namics is far from trivial. In most cases it is done in a phenomenological manner 
by introducing various uncontrolled approximations. 

In this work we address the above problem in the context of the dynamics of 
crystal surfaces. The evolution of crystal surfaces below the roughening transition 
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proceeds by the motion of discrete atomic steps which are separated by high sym- 
metry orientation terraces. One can model step motion by solving the diffusion 
problem of adatoms on the terraces with appropriate boundary conditions at step 
edges. This approach was introduced long ago by Burton, Cabrera and Frank pQ, 
and was further developed by other authors [2]. The resulting models specify the 
normal velocity of each step in the system as a function of its position and shape 
and as a function of position and shapes of neighboring steps. These step flow 
models are capable of describing surface evolution on the mesoscopic scale with 
significant success However, step flow models pose a serious challenge for 

numerical computations, and can be solved only for small systems. 

Several attempts were made to construct continuum models for stepped sur- 
faces [HI El EI El El CHI dl H21 ESI HH CHI CSl HZl CBI , in order to understand their 
large scale properties. The general idea behind these attempts is that step flow 
can be treated continuously in regions where every morphological surface feature 
is composed of many steps. If we label surface steps by the index n, the continuum 
limit in these models is obtained by taking n to be continuous. In what follows we 
will refer to these models as the conventional approach. 

In the literature, there are two methods to derive conventional continuum 
models. One method is to write down the discrete step equations of motion and 
then transform them into a partial differential equation by taking the step index n 
to be continuous |71 ll2lH51ll4U15| . The second method is to start with a continuous 
surface free energy density and derive a surface dynamic equation that minimizes 
it|rfl I%1 ITTI llfil 1171 118j . These two methods are complementary provided that: 1. 
The free energy density of the second method is the continuum analog of the 
free energy of an array of discrete steps. 2. The two methods uses the same mass 
transport mechanism. Such continuum models are fairly successful in describing 
the evolution of smooth surfaces with very simple morphologies. However, they 
suffer from fundamental drawbacks, which do not allow generalizations to more 
complex and realistic situations. 

The most severe drawback is that below the roughening temperature, crys- 
tal surfaces have singularities in the form of corners and macroscopic facets. The 
latter are a manifestation of the cusp singularity of the surface free energy at high 
symmetry crystal orientations. The assumption that every surface feature is com- 
posed of many steps clearly breaks down on macroscopic facets where there are 
no steps at all. Thus, existing continuum models fail conceptually near singular 
regions. Several authors have tried to overcome this problem by solving a contin- 
uum model only in the non-singular parts of the surface and then carefully match 
the boundary conditions at the singular points or lines ^] El El El- In most 
cases however it is not at all clear how these matching conditions can be derived. 
This difficulty is a fundamental drawback of conventional continuum models and 
not merely a technicality. As we argue below, boundary conditions at the singular 
points or lines cannot be derived in the context of conventional continuum models. 

To see why this is true consider the situation near a facet edge. The step at 
the facet edge is special and obeys a unique equation of motion. In contrast to steps 
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in the sloping parts of the surface which all have two neighboring steps of the same 
sign, a facet step has only one neighbor of the same sign and potentially a second 
neighbor of an opposite sign. There might also be special physical conditions such 
as surface reconstruction that add to the uniqueness of a facet step. As we found 
in several cases ^3 ^] 115] . the unique behavior of a facet step sensitively 
determines the amount of material emitted or absorbed at the facet and the rate 
at which steps cross the facet and annihilate. When going to the continuum limit 
these quantities serve as flux boundary conditions at the singularity. The problem 
is that conventional continuum models are derived from the equations of motion 
(or from the surface free energy density) that apply away from the facet and 
are therefore ignorant of the special behavior of facet steps. Thus, the boundary 
conditions at the facet edge must be derived from a careful analysis of the discrete 
dynamics of faces steps. However, in going to the continuum in the conventional 
way, one loses the information regarding the position of individual steps and the 
discrete analysis cannot be performed. 

Another approach for dealing with surface singularities is to round the sur- 
face free energy cusp 1171 H%] , approximating true facets by relatively flat but 
analytic regions. This method avoids the need of specifying explicit boundary 
conditions at the singularity by assuming analyticity of the surface. The correct 
surface behavior is then expected to be captured in the limit of vanishing cusp 
rounding. This procedure completely ignores the key role of facet steps and im- 
plicitly assumes that the surface free energy derived for non singular orientations 
determines the dynamics on facets as well. This assumption is generally false due 
to the same reasons discussed above. An example for a case where cusp rounding 
leads to erroneous results can be found in Appendix 1X1 

In this work we propose a conceptually new definition of the continuum limit, 
which we term Configurational Continuum 19 . Configurational Continuum allows 
construction of continuum models, which arc free of all the limitations of conven- 
tional continuum models discussed above. It provides a rigorous way of deriving 
the continuum model directly from the discrete step equations of motion. Like 
other continuum models, Configurational Continuum has a clear computational 
advantage over the discrete step model due to the small number of discretiza- 
tion points it requires for the description of smooth surface regions in a numerical 
scheme. 



2. Configurational continuum 

In order to overcome the limitations of conventional continuum models we now 
propose a conceptually new definition of the continuum limit for stepped surfaces. 
Our key observation is that a continuous surface height profile can be represented 
by many similar, but not identical, step configurations. Figure ^ is a one dimen- 
sional demonstration of this point. It shows a continuous height function, h{x) 1 
of position x (thick solid line) , and three valid microscopic representations of this 
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profile as step configurations. The main idea of this work is to define the height 
profile in the continuum limit as the upper envelope of the discrete height functions 
of an ensemble of many such step configurations. 




Figure 1 . A schematic illustration of the ensemble of configura- 
tions whose upper envelope defines the continuum limit. 

To construct the configurations of the ensemble let a be the height of a 
single step and N the number of configurations in the ensemble. We construct the 
ensemble so that the height difference Sh, between two adjacent configurations is 
a/N, as depicted in Fig. ^ The continuum limit is obtained when N — > oo. The 
generalization to higher dimensions is straightforward. 

The dynamics of the continuum model is as follows. Each step configuration 
of the ensemble evolves according to the microscopic dynamics. As a result, the 
envelope of discrete height functions changes with time, thus defining the evolution 
of the continuous height function h{r) where r is a vector in the high symmetry 
xy plane. There is one technical complication which might arise if two steps from 
different configurations cross each other. Such an event would make h(r) a multi 
valued function of position and requires a more general mathematical description 
of the surface. For simplicity we ignore this and assume that h(r) remains single 
valued. 

There is a crucial assumption hidden in this definition of the continuum 
model. We postulate that our construction leads to a mathematically well defined 
height function at all times. When does this assumptions hold? Consider two 
initially similar configurations of the ensemble. Our continuum limit is well defined 
provided these two configurations have similar microscopic dynamics and hence 
remain similar at later times. Note however that this assumption has to hold in the 
conventional continuum definition as well, and therefore does not put additional 
restrictions on our model. In fact, if two initially similar configurations evolve 
very differently, one must abandon the continuum limit and follow the specific 
microscopic configuration of interest, using discrete dynamics. 
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We now derive the evolution equation for the continuous height, h(r,t), at 
position r and time t. As a basis for the derivation we assume knowledge of the 
discrete equations of motion for the underlying step flow model. These equations 
of motion specify the normal velocity of a step that passes through (r, f) as a 
function of the local step configuration. We denote this dependency by writing 

v(f,t)=v(C Pit ), (2.1) 

where Cp t t denotes the configuration of steps in the region that influences the step 
velocity at (r , t) . In most models of step flow this region of influence covers a 
small number of neighboring steps. Note that in the context of the discrete step 
model, Cft is the actual configuration of steps in the system. When going to the 
continuum we will be interested in the ensemble of configurations {Cp it } which are 
consistent with h(f,t). 

The continuous height h(r, t) changes with time due to the flow of steps 
through r, and due to nucleation and annihilation of steps. At this stage we dis- 
regard nucleation processes and include them later. First, we consider positions 
which are not local extrema of the height profile. It is obvious from the construc- 
tion of the Configurational Continuum, that for each point r there is exactly one 
configuration t i n the ensemble, which has a step that passes through r at time 
t. That step lies along the unique equal-height contour line, which passes through 
r. As is demonstrated in Fig. ^ the exact positions of neighboring steps in the 
configuration C? tt can be calculated from the knowledge of h(f,t), and the fact 
that in this configuration there is a step at r. Hence, we can use the discrete step 
model Eq. (|2.1(l and calculate the normal velocity of the step v(CV.t). Note that 
at different positions, v (CV,t) is the normal velocity of steps which may belong to 
different configurations in the ensemble. 

Next we define the directional gradient in the direction from which steps flow 
towards r 

s -/^ % h(r — e-v(r,t)) — h(r,t) ,„ „. 

Vh-df,t) = -v(r,t) lim — i , (2.2) 

where v(r,t) = '*), . This is useful for the calculation of the current of steps 
arriving at f: 

N 

J{r,t) = - V/i_«(f,t) ■ v(C f , t ) , (2.3) 
a 

where we have used the fact that the local step density is |V/i_«(r, t)\ N/a. Note 
that J is the current of steps belonging to all configurations in the ensemble, 
and not to one particular configuration. Since each step (from any configuration), 
which passes through f changes the height of the ensemble envelope by a/N, the 
continuous height profile obeys the evolution equation 
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The above derivation of the evolution in the continuum is not valid at local 
extrema of the surface, because generally one cannot define a unique equal-height 
contour line which passes through such a point. To avoid the problem, we define 
dh/dt at local extrema as the limit of the height time derivative as one approaches 
these points. This limiting procedure is justified, since there are no microscopic 
realizations of the surface with steps exactly at the local extrema. 

At this point we emphasize that the Configurational Continuum evolution is 
formally identical to the evolution of the discrete step model. This statement is 
almost trivial, since the definition of Eq. I|2.4(l follows the envelope of the ensemble 
of configurations, and each configuration evolves with step velocities calculated 
from the discrete step model. Thus Eq. 12.41) is exact. Moreover, the use of di- 
rectional derivatives in the derivation of Eq. I|2.4fl makes it valid even at singular 
surface regions such as corners or facet edges. Similarly to other continuum models 
it is solved numerically by discretization of space, which is the only approximation 
involved in such solutions. 

What is the relation between Configurational Continuum and conventional 
continuum models? In regions where h(f, t) is analytic, the evolution equation (|2.4|) 
reduces to the continuity equation 



If h(r, t) is sufficiently smooth, v (Cj?t) can be approximated as a local function of 
h and its spatial derivatives, as is commonly done in conventional continuum mod- 
els. Making this approximation will therefore recover the conventional continuum 
approach. We can conclude that in analytic surface regions the conventional con- 
tinuum approach approximates the Configurational Continuum model and that 
the approximation quality depends on the smoothness of the surface. However, 
near corners, facets or regions where the profile is not smooth, one cannot recon- 
struct the microscopic step configuration from the local value of h and its spatial 
derivatives. In these regions v (Cp.t) contains non local information and as a result 
Eq. H2.4fl cannot even be written as a differential equation. 

Is there any computational gain in using such a continuum model? After 
all, we replaced a discrete model, which follows the evolution of a single micro- 
scopic step configuration, by a model which follows a whole ensemble of step 
configurations. The key point is that we do not have to follow all the steps of all 
configurations. To calculate dh(r,t)/dt, it is enough to locally follow the single 
configuration which has a step that passes trough f at time t. In addition, the 
continuum evolution equation is solved on a grid, and the density of grid points 
can be very small in regions where the continuous height profile is smooth. The 
smoothness of the profile allows very accurate interpolation between these points. 
Only near singular points or lines we have to use a rather dense grid, and there 
is no computational gain in these regions. In practice, the total number of grid 
points used can be orders of magnitude smaller than the number of points one has 
to use in order to follow the evolution of a single microscopic step configuration. 



dh(r, t) 
dt 



(2.5) 
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So far wc ignored the possibility of island or void nucleation. It is possible to 
include island or void nucleation in our model provided that we have a microscopic 
description for these events which determines the nucleation probability in a given 
step configuration. Within our continuum approach, the nucleation probability 
at a point on the continuous surface is the ensemble average of the microscopic 
nucleation probabilities at this point. For demonstration proposes we consider a 
simple model where the probability for the nucleation of an island on a terrace 
grows as the square of the local concentration of diffusing adatoms. Information 
regarding the values of terrace adatom concentrations is already contained in the 
underlying step flow model Eq. JOJ, since it is used in the calculation of adatom 
fluxes into and out of steps. 



3. Numerical solutions of the Configurational Continuum 

We now apply our approach to a few simple cases. First, we consider a conic 
structure which consists of circular concentric steps. This crystalline cone was 
studied in Rcfs. |12U13) where we wrote a one-dimensional step flow model for the 
step radii in the absence of growth flux. In the diffusion limited case where adatom 
diffusion across terraces is the rate limiting process, the equation of motion for the 
step radii r n read: 



dt ksTri In (n/ra) 

dr n D S C^Q, ( \i n - fi n+ i fin-i - fi n 



dt k B Tr n \\n(r n /r n+ i) hx(r n -i/r n ) 

with the step chemical potentials /i„ given by 

2r„+i 1 2r n _i 



nr 

[i n = h 

r„ 



r n +i+r n {r n + x -r n f r n + r n _i (r n - r n ^) 3 



n>l, (3.1) 



(3.2) 



In the above expressions D s is the adatom diffusion constant, f2 is the atomic area 
of the crystal and is the equilibrium concentration in the vicinity of a straight 
isolated step. T is the step line tension, G is the step-step interaction strength, 
T is the temperature and Ub is the Boltzmann constant. Eq. (|3.2I) can be used 
for calculating the chemical potential of the top step, fj,i, by omitting the second 
interaction term. 

Numerical integration of the above step model shows that the cone decays 
through the collapse of individual islands. During the decay a facet develops at 
the top of the cone with a radius that grows with time as t 1 / 4 . Similar equations 
can be written in the presence of growth flux. With flux, the cone grows except 
at the peak which initially decays and then saturates. A facet forms at the peak 
after saturation. 

The continuum model we solved in the context of this example is a fully 
two dimensional model, which can, in principle, develop non radially symmetric 
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morphologies. The microscopic dynamics we used were a two dimensional general- 
ization of the microscopic equations for step radii of the discrete one dimensional 
model given above. In particular the step chemical potential was generalized to 



where n (r) is the local step curvature. r<j and r u are the coordinates at which the 
lower and upper neighboring steps are closest to the step at r. It is assumed that 
these neighbors have the same sign as the step at r and that steps of opposite 
signs do not interact (in which case the relevant interaction terms in Eq. (|3.3|1 
are omitted). This expression was derived assuming an l~ 2 repulsion between two 
segments of two different steps which are separated by a distance I. Under this 
assumption, Eq. I|3.3|l is exact to first order in the curvature of neighboring steps 
and qualitatively captures the interaction when they have large curvatures. 

The equations describing adatom diffusion on terraces were solved assuming 
the steps at f, and f u are circular and concentric and that the radius of the 
step at r is l/fc(r). Any microscopic dynamics, such as a full solution of the dif- 
fusion equation on each terrace, or a more detailed calculation of the interactions 
between steps can easily be used in the framework of our model. For the sake of 
demonstrating the validity of Configurational Continuum the simple dynamics we 
chose are sufficient. 

Figure [21 shows a comparison between a numerical solution of the discrete 
step model and a cross section from the two dimensional solution of our contin- 
uum model in the absence of growth flux and when nucleation of new steps is not 
allowed. Clearly the continuum model captures the main features of the surface 
evolution. In particular the width and height of the top facet are in excellent agree- 
ment with the discrete model. Fig.[3]shows a similar comparison in the presence of 
growth flux. Again the agreement is quite impressive. In this case we allowed new 
islands to nucleate and our simulation indicates that nucleation events occur on 
the top facet once it becomes large enough. There is hardly any nucleation on the 
finite slope regions, because the steps there efficiently absorb the excess material. 
Figure|3Jc) shows an island on the top facet, which nucleated and started growing. 

We now turn to the more demanding example of bidirectional sinusoidal 
grating relaxation. Here the initial surface height profile of wave length L is given 



The relaxation of this profile towards a flat surface was studied by Rettori and 
Villain |2(J| . who gave an approximate solution to a step flow model, in the limit 
where the interaction between steps can be neglected with respect to the step line 
tension. The decay of bidirectional sinusoidal profiles was also studied numerically 




(3.3) 



by 
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Figure 2. Decay of a crystalline cone. The solid line shows a 
cross section of the two dimensional solution of Eq. (|2.4(l . Circles 
show the surface evolution according to the one dimensional step 
flow model Eq. (|A.9|I . (a) is the initial morphology and (b) and 
(c) show the surface at later times. 



|18j . We now apply our model to this problem assuming diffusion limited kinetics 
without deposition flux or island nucleation. 

For weak interactions between steps, the surface height evolves according to 
dh^/dt oc V 2 k ~ L~ 3 , where n is the step curvature. We therefore expect the 
following scaling law for hi,(f, t): 

h L (r,t) = h L=1 (r/L,t/L 3 ) . (3.4) 

Figure 0] shows the data collapse of cross sections of profiles resulting from 
our continuum model. The different symbols correspond to profiles of different 
wave lengths at time t = toL 3 for some fixed to. The quality of the data collapse 
shows that the scaling scenario (|3.4II holds very accurately. Note that large facets 
have developed at the surface extrema, and they are connected by very steep 
slopes. This shape does not agree with Rettori and Villain's heuristic argument 
|2(J|. which predicts that facets appear also near h = lines. Nevertheless, their 
prediction that after a short transient the amplitude of the height profile decays 
as t/L 3 is in agreement with both Eq. I|3.4[l and our numerical solutions. 

Figure El shows results for a much stronger repulsive interactions between 
steps, where the scaling law (|3.4f) clearly does not hold. Fig. [3] (a) shows profiles 
of different wavelengths which have relaxed to half of the initial amplitude. Here 
profiles with a smaller wavelength have small facets and moderate slopes. This 
happens because repulsion between steps becomes increasingly important as the 
profile wavelength is reduced. At long wavelengths the weak step-step interaction 
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Figure 3. Evolution ol a crystalline cone under growth condi- 
tions. The lines show cross sections of the two dimensional solu- 
tion of Eq. (|2.4() . Circles show the surface evolution according to 
the one dimensional step flow model. The initial shape is the ex- 
act cone shown in Fig. [21 (a), (a) and (b) show the generation and 
evolution of a facet at the top at later times, (c) shows an island 
which nucleated (solid) and grew (dashed) on the top facet. Here 
no comparison is made with the step flow model. 

limit of Fig. 0| is recovered. Fig. 03(b) shows the different amplitudes as a function 
of scaled time t/L 3 . 

4. Conclusions 

We proposed a novel continuum model for the evolution of sub roughening crystal 
surfaces. Our model, which we term Configurational Continuum, is derived directly 
from the underlying dynamics of atomic steps. Unlike conventional continuum 
models, Configurational Continuum is fully consistent with the step dynamics and 
is capable of handling singular surface features such as facets and corners. 

The key idea in our model is to view the continuous surface profile as the 
envelope of an ensemble of step configurations which are all consistent with the 
continuous profile. Knowing the ensemble envelope, it is always possible to re- 
construct individual configurations and evolve them in time. This evolution of 
individual configurations determines the evolution of the ensemble envelope which 
within our model is interpreted as the evolution of the profile. The continuum 
limit in our model is naturally realized because the continuous profile induces a 
continuum of possible step configurations. 
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Figure 4. Data collapse in the evolution of bidirectional sinu- 
soidal profiles with different wave lengths when the repulsive in- 
teraction between steps is very weak. Wave lengths shown are 
L = 64 (circles), L = 128 (squares) and L — 256 (triangles), a) 
Cross sections of the profiles with different wave lengths measured 
at time t = t n L 3 . The cross sections are along the y = —L/A line 
(peak to valley) . b) Amplitude decay of the different wave lengths 
as a function of scaled time t/L 3 . 



Like other continuum models, Configurational Continuum has a computa- 
tional advantage over the underlying discrete step model. When solved on a com- 
puter, it is possible to use a sparse grid in regions where the profile is very smooth. 
However, near corners or facets our model requires a fine discretization grid. The 
fine grid is necessary in order to faithfully reconstruct the step configurations which 
are possible near a singular point. Our model thus has the important property of 
being capable of describing step flow on different scales in a consistent way. In 
smooth surface regions, Configurational Continuum provides a coarse grained nu- 
merical description of surface dynamics. However it is still capable of accounting 
for the unique behavior of steps near singular points or lines. 

The problem of connecting between different scales in dynamical systems is 
not limited to the evolution of surfaces. This problem is widespread in physics, 
engineering and biology as well as in other fields. Our hope is that the ideas at the 
basis of Configurational Continuum can be applied in other multi-scale problems. 
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Figure 5. Bidirectional sinusoidal profiles in the case of stronger 
step-step interactions. Wave lengths shown are L — 64 (solid), 
L = 128 (dashed) and L — 192 (dotted), a) Cross sections of the 
profiles with different wave lengths which have relaxed to half of 
their initial amplitude. The cross sections are along the y = —L/A 
line (peak to valley), b) Amplitude decay of the different wave 
lengths as a function of scaled time t/L 3 . 



Appendix A. Failure of the cusp rounding method. 

In this appendix we give an example which shows how rounding of the surface free 
energy cusp can lead to erroneous solutions for surface evolution. As an example 
system we again choose the crystalline cone studied in Refs. ^JEj- For simplicity 
we consider the diffusion limited case. This system is convenient since it exhibits 
scaling. In the scaling limit conventional continuum modelling of the cone becomes 
exact |12l I13| and we can concentrate on effects introduces by the cusp rounding 
method. We start by solving a model with a rounded free energy cusp and later 
compare the solution of this model with the relevant discrete step model and with 
the solution of the Configurational Continuum. 

In the conventional continuum approach, the continuous free energy density 
of sub roughening surfaces has a cusp singularity at the high symmetry orientation. 
In a coordinate system where (x, y) is the high symmetry plane and h(x, y) is the 
surface profile, the projected (on (x,y)) surface free energy density assumes the 
form |21j : 

T{x,y) = T Q + r\Vh(x iy )\ + j\Vh(x,y)\ 3 . (A.l) 
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This form is the continuum analog of the free energy of an array of steps with line 
tension T and an inverse square repulsive step-step interaction of strength G. 

The singular nature of the surface free energy complicates the modelling of 
surface evolution. Several authors have tried to overcome this problem by rounding 
the cusp with a small parameter a: 



T a (x,y) = T Q + T \{yh{x,y)f + c?\ +^[(Wh(x,y)Y + a A \' . (A.2) 

The hope behind this regularization scheme is that in the limit a — > the resulting 
model captures the correct surface dynamics. 

Surface dynamics is derived from T a as follows. Taking the functional deriv- 
ative of J- a we obtain the surface chemical potential: 

Sh \dx dh x dy dh y J 

where h x — dh/dx and h y = dh/dy. For a radially symmetric profile h(r,t) we 
find that 

n 



1/2 G 



\ 2 , „2 



3/2 



V + h rr) (r + G {hi + « 2 )) + h% r (G j^—^ 



(A.4) 

where h r = dh/dr and h rr = d 2 h/dr 2 . 

In diffusion limited kinetics variations in the chemical potential give rise to 
currents which are proportional to the chemical potential gradient. For our radially 
symmetric profile we can consider only the radial component of this current 

The dynamic equation for the profile can now be written using the continuity 
equation 

dh np s c^ i d ( d^ a \ 

— = -n ■ VJ = - - -— r— . (A.6) 



dt ksT r dr \ dr 

The a — > limit of Eq. IjA.Bj) automatically gives the correct equilibrium 
crystal shapes because Eq. (|A.2J| goes to Eq. (|A.1|) in this limit. We want to check 
whether this limit also gives the correct (consistent with step flow) dynamics. 

We applied Eq. (|A.6|I to a crystalline cone h(r, t = 0) ~ — r. The similarity 
sign here indicates that the tip of the initial profile was smoothed in order to have 
an analytic surface. Analyticity at the origin was also used as a boundary condition 
for the surface evolution. Numerical solutions show that, at long times the profile 
slope obeys a scaling law 

h r (r,t) ^-Firr 1 ^) . (A.7) 

This behavior agrees with the scaling properties exhibited by a discrete step flow 
model of the same surface structure ^JE]- In Fig. El we show the resulting scaled 
slopes (dots) for different values of the cusp rounding parameter a. These long time 
solutions are not sensitive to the initial smoothing of the cone. As a is reduced 
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we observe the appearance of a very fiat region around the origin. This flat region 
supposedly becomes a true facet in the a = limit. 

For small values of a, solutions of the dynamic equation (|A.6(I approach 
scaling very slowly. For this reason it becomes increasingly difficult to probe the 
a = scaling state. In order to reach smaller values of a we continued in the 
following way. We assumed that the scaling ansatz l|A.7(l holds and used it to 
transform Eq. (|A.6() into an ordinary differential equation for the scaling function 
F(£) in the scaling variable £ = ri -1 / 4 . Replacing h r (r,t) in Eq. (|A.6|I by the 
scaling function we obtain the following equation: 



1 



F' 



VtD s C eq d 
k B T ~dj 
Q 



1 d 
Id! 
F 

r 

+F 2 F' 



d£ 

F'^j (r + G (F 2 
G r 



2 )) 



F 2 



(A. 



with F' = dF/d£. 

Solutions of the above equation for the large values of a agree with the scaling 
states of the dynamic equation l|A.6|) . The dashed lines in Fig.Elshow the resulting 
scaling functions for smaller values of a. Finally in order to determine the a = 
limit we estimated the scaled position of the facet edge from the a^O solutions. 
This position selects [T2l IT3] a unique scaling solution for the a = case of Eq. 
(|A.8J) . The resulting scaling function is shown by the solid line in Fig. By our 
procedure this function is an approximation for the true scaling function of the 
system according to the cusp rounding method. 

The solid line in Fig. should be compared with the behavior of a system of 
discrete steps. For this purpose we introduce the following step model: 
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This is the same step model studied in section|3and in Refs. |12U18j with modified 
step chemical potentials. The modification has a small effect on the model behavior 
and does not introduce any qualitative changes. In particular, this step model obeys 
the same scaling properties that were studied in Refs. |12M13| . i.e., the density of 
steps in this model scales according to D(r,t) — F^ iscrete {rt _1 / 4 ). In addition, 
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Figure 6. Scaling functions of the slope for different values of the 
cusp rounding parameter. Dots show scaling functions obtained 
from a direct solution of the dynamic equation (| A. 6|> . Dashed lines 
show solutions of the scaling equation (|A.8|I . The solid line shows 
the estimated a = solution of Eq. (|A.8Jl . 



applying the scaling analysis of Refs. ^JEH to this modified step model results in 
an ordinary differential equation for the scaling function Fdiscrete which is identical 
to the a = limit of Eq. {Qj) . This fact gives us a basis for comparison between 
the step model and the cusp rounding scheme. Identifying the step density of the 
discrete model with the slope of the continuum model we can finally compare the 
limiting solution from the cusp rounding method with the scaled density of steps. 
In Fig. [7| we show that these two functions do not coincide. The scaled position 
of the facet edge in the cusp rounding method is about 40% larger than the one 
in the discrete system. This means that the difference between the two facets in 
real space diverges at long times. The height at the origin according to the cusp 
rounding method will suffer from the same errors. By assuming analyticity of the 
profile throughout the limiting procedure of the cusp rounding method, we have 
imposed erroneous boundary conditions at the facet edge. 

Figure[7|also shows a one dimensional solution of the Configurational Contin- 
uum model for this cone system. The Configurational Continuum model predicts 
scaling of the slope as well. The discrepancy between the resulting scaling func- 
tion and the discrete step system is much smaller and is consistent with what one 
would expect from discretization errors. 
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Figure 7. Comparison between the scaling function predicted by 
the cusp rounding procedure (solid), scaling function from a one 
dimensional solution of the Configurational Continuum (dashed) 
and the scaled density of steps in the discrete model (circles). 
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